Bound-state energy of the d = 3 Ising model in the 
broken-symmetry phase: Suppressed finite-size corrections 

Yoshihiro Nishiyama 
Department of Physics, Faculty of Science, 
Okayama University, Okayama 700-8530, Japan 
(Dated: April 10, 2008) 

Abstract 

The low-lying spectrum of the three-dimensional Ising model is investigated numerically; we 
made use of an equivalence between the excitation gap and the reciprocal correlation length. In 
the broken-symmetry phase, the magnetic excitations are attractive, forming a bound state with an 
excitation gap m,2{< 2m\) (mf. elementary excitation gap). It is expected that the ratio m^jfnx is 
a universal constant in the vicinity of the critical point. In order to estimate m^/mi, we perform 
the numerical diagonalization for finite clusters with N < 15 spins. In order to reduce the finite- 
size errors, we incorporated the extended (next-nearest-neighbor and four-spin) interactions. As a 
result, we estimate the mass-gap ratio as milm\ = 1.84(3). 



I. INTRODUCTION 



The first-excitation gap (reciprocal correlation length) is of fundamental significance in 
the theory of critical phenomena. On one hand, little attention has been paid to the low-lying 
spectrum. Recent consideration on the Ising model, however, revealed rich characters in the 
low-lying spectrum . Actually, in the broken-symmetry phase, 

the magnetic excitations are attractive, forming a pair of bound states with characteristic 
excitation gap {rrii}. 

For example, in two dimensions (d = 2), there appear eight types of excitations (i = 
1,2, ... ,8) l|, [2|. Each of these excitations possesses an intrinsic mass gap m.j. That is, 
for each excitation, the ratio rrtij 'mi (mj: an elementary excitation gap) is a universal 
constant in the vicinity of the critical point. [Notably enough, the elementary excitation 
m 1 is also a composite particle (bound state) |3j. In this sense, all particles are equally 
elementary, reflecting a highly non-perturbative nature of this problem.] The point is that 
in the continuum (scaling) limit, the d = 2 Ising model is exactly solvable even in the 
presence of the magnetic field; see Ref. [3j] for a review. As mentioned afterward, (properly 
scaled) infinitesimal magnetic field is significant to stabilize the bound state. 

On the contrary, in d = 3, such an exact solution is not available. The detailed structure 
of the d = 3 spectrum is not fully understood. So far^the lowest bound state m 2 , namely, 

3, H, S, H, Q, 



the second-excitation gap, has 



Deen studied in detail 



states is considered in Refs. , [lO, For example, by means of the Monte Carlo 

method [4], a clear indication of a bound state m 2 was found. The simulation data indicate 
the mass-gap ratio 

ma/mi = 1.83(3). (1) 

(See also Refs. .) Note that the Monte Carlo method does not yield the excitation 

gap m 2 directly. Rather, one has to look into the asymptotic form of the correlation function, 
and estimate the sub-dominant decay rate oc l/m 2 together with the dominant one cx \jm\. 
(Note that the decay rate is related to the inverse of the corresponding mass gap.) On the 
contrary, the numerical diagonalization method is suitable for investigating the spectrum 
{^j). However, in practice, the tractable system size with the diagonalization method is 
severely restricted; the numerical diagonalization requires a huge computer-memory space. 
To cope with this difficulty, in Ref. [9( , the authors truncated the number of Hilbert-space 



lOj ; a series of bound 



2 



bases (stochastic diagonalization method) to calculate the spectrum approximately. The 
result appears to be consistent with the above-mentioned Monte Carlo result. (However, 

n 

as indicated in Table 1 of Ref. |9j], the results m 2 /m 1 scatter, depending on the distance 
from the critical point. In this paper, we carry out systematic finite-size-scaling analysis to 
resolve such uncertainty.) 

In this paper, we employ the numerical diagonalization (transfer matrix) method to 
calculate 11x2! vn\. In order to cope with the difficulties mentioned above, we make the 
following modifications. First, we improve the finite-size-scaling behavior. That is, we 
reduce corrections to finite-size scaling by tuning plaquette-type (next-nearest-neighbor and 
four-spin) interactions. Such a reduction of corrections for the d = 3 Ising model was pursued 



in Refs. 



12,ll3|,ll4(, where the criticality of the d = 3 Ising model was investigated. (Similar 



attempt was made in a lattice- field-theoretica 
we follow the formalism advocated in Ref. 



context 



16 



13,[l8|,[l9j.) To be specific, 



141 ]. where the Hamiltonian reads 



H — —Jnn SjSj — Jnnn SjSj — Jo SiSjSkSi — Si. (2) 

(i,J> «M» [i,3,k,t\ i 

Here, the Ising spins {Si = ±1} are placed at the cubic-lattice points i, and the summations, 
S/i j\i Ylm 3))! an d j,kih run over an nearest-neighbor pairs, next-nearest-neighbor spins, 
and round-a-plaquette spins, respectively. The parameters {J a } are the corresponding cou- 
pling constants. As usual, H denotes the uniform magnetic field. [We set the temperature 
as a unit of energy (T = 1) throughout this paper.] We survey the one-parameter subspace 

(Jnn, Jnnn, Jo) — Jnn 1 1, , -~ — J , (3) 

\ Jnn Jnn J 

containing a renormalization-group fixed point, 

(Jnn, Jnnn, Jn) = (0.1089828666435, 0.0445777727956, -0.0065117950492), (4) 

at Jnn = Jnn', the fixed point was determined approximately for a real-space-decimation 
transformation 114]. Around the fixed point, the finite-size corrections may cancel out be- 

nn 

cause of the absence of irrelevant interactions [13|, HJ]. Second, we carry out two types of 
finite-size-scaling analyses for H 7^ and H = 0. Such independent information allows 
us to estimate m^/mi systematically. Last, we utilize Novotny's transfer-matrix method 



2d, 



21], which enables us to treat an arbitrary (integral) number of spins iV = 5, 6, . . . , 15 



3 



constituting a unit of the transfer-matrix slice; note that conventionally, the number of spins 
is restricted within N — 4, 9, 16, . . . . 



The ratio m 2 /mi was considered analytically in Ref. 1(|. In the following, we outline the 
argument in order to elucidate the underlying physics. The authors calculate the four-point 
vertex (Bethe-Salpeter equation) of the 4 model perturbatively up to the second order; the 
analyticity of vertex contains information on the bound state. As a result, they expressed 
1712/1711 in terms of the coupling constant at the critical point; the critical coupling itself 
is a nontrivial parameter. By using the critical coupling 22] determined numerically, the 
authors arrive at a convincing result 1122/1711 = 1.828(3) at the leading order. However, at 
the next-leading order, the ratio m 2 /m! turned out to be negative, indicating that m 2 /m! 
is highly non-perturbative. Nevertheless, their consideration reveals a close relationship 
between 7712/ vn\ and the critical coupling (critical amplitude). 

In fairness, it has to be mentioned that the spectral structure of the d = 3 Ising model 



is considered in Refs. 



23 



24 



251 ] . Here, the relationship between the scaling dimension 



and the corresponding excitation gap is investigated for various geometries and boundary 
conditions via the conformal invariance. In the present paper, we dwell on the off-critical 
excitations having a physical interpretation as magnons. 

The rest of this paper is organized as follows. In Sec. [Ill we present the simulation 
results for m 2 /m 1 . We also explain the simulation scheme. In the last section, we present 
the summary and discussions. 



II. 



NUMERICAL RESULTS 



In this section, we present the numerical results for the d = 3 Ising model, Eq. ([2]), w ith 
the extended interactions, Eq. ([3]). We employ Novotny's transfer-matrix method 20l. 121] 



for finite clusters with N = 5, 6, . . . , 15 spins. To begin with, we outline the simulation 
algorithm. 



A. Simulation method: Novotny's technique 



In this section, we explain the simulation scheme 



consult with Ref. 



2l) briefly; we refer readers to 



211 ] for details. 
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A schematic drawing of a unit of the transfer-matrix slice (cross-section of the transfer- 
matrix bar) is presented in Fig. [TJ Basically, the spins constitute a d — 1-dimensional 
(zigzag) alignment {Si} (i = 1,2, . . . , N). The dimensionality is lifted to d = 2 by the 
bridges over the (iV 1//2 )th neighbor pairs. Because the basic structure is one-dimensional, we 
are able to construct the transfer-matrix unit with an arbitrary (integral) number of spins; 
note that conventionally, the number of spins for the d = 2 structure is restricted within 
N = 4, 9, 16, . . . (quadratic numbers). The linear dimension L is given by 

L = N l '\ (5) 

because the N spins form a rectangular cluster. We treated the system sizes N = 5, 6, . . . , 15 
with this technique. 

Our main concern is to calculate the excitation gap. We perform the numerical diago- 
nalization within the subspace k = (k: the wave number). Additionally, we consider the 
spin-inversion-symmetry index ± in the case of H = 0. First, for H ^ 0, the zth excitation 
gap is given by the formula 

rrii = ln(A /Ai). (6) 

Here, the parameters {Aj} are the eigenvalues of the transfer matrix, satisfying the inequality, 
Ao > Ai > A2 > • • • • Second, in the case of H = 0, the first and second excitation gaps are 
given by 

mi = hi(A +/Ar) (7) 

and 

m 2 = ln(A+/A+), (8) 

respectively, with the transfer- matrix eigenvalues {A^} specified by the spin- inversion- 
symmetry index ±. 

In practice, we need to modify Novotny's method in order to incorporate the extended 



interactions such as the next-nearest-neighbor and four-spin couplings 20j. We refer readers 



to consult with Ref. [2lJ for a detailed description of the algorithm. 

Lastly, to avoid confusion, we address a remark on the fixed point, Eq. (j3J). The fixed 
point is not used to analyze the criticality such as the mass-gap ratio mijm\. We utilize the 
conventional finite-size-scaling method. We perform the computer simulation around the 
fixed point, where corrections to scaling cancel out satisfactorily; see the Introduction. In 
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this sense, the fixed-point analysis is a preliminary to the finite-size scaling. Owing to the 
suppression of finite-size errors, we are able to analyze the simulation data systematically. 



B. Finite-size-scaling analysis of mi/mi'. H ^ 

In this section, we analyze the simulation data under an infinitesimal magnetic field 
H ^ 0. The magnetic field stabilizes the bound state. (The case of H = is considered in 
the next section.) 

To begin with, we consider the finite-size scaling for the mass-gap ratio 

— = / ((Jnn - J* NN )L l '\ HW») , (9) 
rrii 

with the critical point J NN = 0.11059 [jjl]. The critical indices, v and t/h, are taken from 
the Monte Carlo estimates v = 0.63020 and yn = 2.4816 (3d-Ising universality) reported in 



Ref. |26j. The above scaling relation, Eq. (Q, stems from the fact that the mass-gap ratio 
is a dimensionless (scale invariant) quantity. 

Based on the above scaling relation ([9]), in Fig. [2} we presented a plot, ( Jnn — Jnn)^ 1 ^- 
7712/7%, for various Jnn, N = 5, 6,..., 15, and the fixed scaling parameter HL Vh = 4; 
afterward, we argue the constraint, HL Vh = 4, in detail. We see that the data collapse into 
a scaling curve. We stress that there are no fitting parameters in the finite-size scaling of 
Fig. [21 Actually, all scaling parameters, J NN , v, and y^, are taken from the existing values, 
as mentioned above. 

Noticeably enough, there appears a plateau with 1712 /mi 1.8 extending in the low- 
temperature side, Jnn — Jnn > 0- This feature indicates an existence of a bound state with 
777.2/7711 ~ 1.8 in the broken-symmetry phase. As mentioned in the Introduction, a properly 
scaled magnetic field is useful to observe the bound state clearly. 

The plateau of Fig. [2]yields the mass-gap ratio. We read off m 2 /mi at (Jnn — Jnn)L 1 ' v = 
0.3 for each system size, and plotted the values against 1/L 2 in Fig. [3l (Afterward, we 
consider the validity of the abscissa scale 1/L 2 .) The least-squares fit to the data for N > 7 
yields an estimate m 2 /mi = 1.835(4) in the thermodynamic limit. 

Lastly, we address a few remarks. First, we consider the scaled magnetic field H = 4/L Vh . 
Owing to the scaling, the second argument of the scaling relation, Eq. (Q, is kept invariant; 
namely, the presence of magnetic field becomes implicit under the scaling. On one hand, the 
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magnetic field plays a significant role to stabilize the bound state. Second, we explain why 
the coefficient A = 4 in the scaling H = A/ L Vh is an optimal one. We surveyed various values 
of A, and found that for A > 5, the plateau inclines; namely, it becomes ambiguous where 
we measure the plateau hight. For A < 3, eventually, there occur successive level crossings, 
smearing out the plateau structure; the magnetic field H may be too weak to stabilize the 
plateau (bound state). Last, we consider the extrapolation scheme 1/L 2 (abscissa scale in 
Fig. [3]). As shown in the figure, the data for H ^ 0, and H = are consistent with 
each other. Moreover, using the v data reported in Ref. [fjj], we arrive at a convincing 
result v = 0.6314(14) through the 1/L 2 extrapolation scheme; the result is comparable with 
that of Monte Carlo v = 0.63020(12) 26j. (On the contrary, in Ref. 14|, a rather naive 
extrapolation scheme 1/L was used, resulting in v = 0.6245(28).) In principle, the finite-size 
corrections should obey the formula 1/L U with uj = 0.821(5) [26J. That is, the convergence of 
our result is rather accelerated. This fact indicates that the finely tuned coupling constants, 
Eq. ([3]), suppress the finite-size corrections significantly. So far, the technique has been 
applied to investigating the equilibrium properties [l3|, ll4j. Our result indicates that the 
technique is applicable to the spectral property as well. 



C. Finite-size-scaling analysis of 7712/ 'mi'. H = 

In this section, we consider the case of H = 0. Here, we have to make the following 
replacement 

m 3 /m 2 — » m 2 /mi, (10) 

because the ground state is doubly degenerated owing to the spontaneous symmetry break- 
ing; namely, the first-excitation gap m\ = does not make any sense. (More specifically, the 
excitations are given by Eqs. (J7|) and (JED-) Actually, to avoid this complication, we applied 
a magnetic field, and resolve the degeneracy in the above section. 

In Fig. HI we present the scaling plot, (Jnn — Jnn)^' 1 ^' 111 ^/ 111 !-, f° r various Jnn, N = 
5, 6, . . . , 15 and the fixed magnetic field H = 0. The scaling parameters are the same as 
those of Fig. [2j The data for > 7 collapse into a scaling curve satisfactorily. 

From Fig. HI we read off m 2 /mi at (Jnn — Jnn)^ 1 = (critical point) for each N, and 
plotted them against 1/L 2 in Fig. [3j The least-squares fit to the data for N > 7 yields an 
estimate 7712/ mi = 1.841(4) in the thermodynamic limit. The result is consistent with that 
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of H ^ (Sec. Ill Bp . In the last section, we address a concluding remark. 



D. Scaling analysis of the ordinary Ising model Jnnn = </□ = 

In this section, tentatively, we turn off the extended interactions, Jnnn — Ju — 0. That 
is, we consider the ordinary Ising model, which is relevant to the preceding Monte Carlo 



study [l]. 

In Fig. El we present the plot, (Jnn — Jnn)^ 1 m 'i/ m i, with J NN = 0.22165455 26], 
N = 5,6,.. .,15, Jnn = Ju = 0, and H = 0; the other scaling parameters are the same as 
those of Fig. HJ 

We see that the data scatter as compared to those of Fig. HJ where we incorporated 
the extended interactions. This result demonstrates that the extended interactions suppress 
the finite-size errors significantly. Such a suppression of corrections was observed in Refs. 







ljj], where the thermal equilibrium properties, rather than the spectral ones, were 
explored with this technique. 

As mentioned above, the Monte Carlo simulation [4] was performed at H = and Jnn = 
Ju — 0. We argue an implication of the present simulation data as to the preceding Monte 
Carlo study. Our scaling analysis, Fig. El indicates that the scaling regime, where we observe 
correct mass-gap ratio mijmx ~ 1.8, is extremely narrow (Jnn — Jnn)^ 1 ^ < 0-2- On the 
other hand, the Monte Carlo simulation was performed at ( Jnn — Jnn)^ 1 ^ ~ 2-16 and 2.46 
with (Jnn, L) = (0.23142, 30) and (0.2275, 45), respectively; those parameters are rather out 
of the scaling regime. Nevertheless, in such close vicinity to the critical point, the notorious 
critical-slowing-down problem would emerge. Possibly, the data under H 7^ may provide 
valuable information on mijmx as demonstrated in Sec. Ill Bl 

III. SUMMARY AND DISCUSSIONS 

The bound-state energy mi of the d = 3 Ising model in the broken-symmetry phase, 
Jnn — Jnn > 0, is investigated numerically. The mass-gap ratio mijmx (mx' elementary 
excitation gap) is a universal constant in the vicinity of the critical point. We employed 



Novotny's transfer-matrix method 



20 



2l| to treat a variety of system sizes iV = 5,6,...,15. 



In order to improve the finite-size-scaling behavior, we extended the interactions, Eq. 
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following Ref. Owing to this improvement, we attain improved finite-size corrections 

(Fig. H]), as compared to those of the ordinary Ising model (Fig. EJ); 

In Figs. [2] and HI we made different approaches to the critical point through H ^ and 
H = 0, respectively. As shown in Fig. [31 these results provide lower and upper bounds for 
7^2 A^i; respectively. Combining these results, we estimate the mass-gap ratio 

m 2 /m 1 = 1.84(3). (11) 



The result is consistent with that of Monte Carlo m-il'mx = 1.83(3) [4]. However, as argued 
in Sec. Ill CI there arises a subtle discrepancy as to the scaling regime. 

So far, the technique to suppress finite-size corrections has been utilized to survey the 
thermal-equilibrium (ground state) properties fl Q . The present analysis demonstrated 
that this technique is also applicable to the excitation spectrum. Likewise, the technique 
would be of use to investigate the dynamical properties such as the spectral intensity and 
the relaxation to thermal equilibrium. This problem is addressed in future study. 
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FIG. 1: A schematic drawing of a transfer-matrix unit of the d = 3 Ising model (|2|). As indicated, 
the spins constitute a d = 1-dimensional (zigzag) alignment {Si} (i = 1,2, ... , N), and the dimen- 
sionality is lifted to d = 2 by the bridges over the (iVV2)th -neighbor pairs. (A transfer-matrix unit 



for d = 3 should have a d = 2 structure.) Full details of the algorithm are presented in Ref. 
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FIG. 2: Finite-size-scaling plot, (Jaw — J^ N )L l l l '-rr^/mi, for the (i = 3 Ising model under the 
scaled magnetic field H = A/L Vh is shown. Here, we set J^ N = 0.11059 [lj], and {v,yh) = 
(0.63020,2.4816) (3d- Ising universality) [26]. From the plateau in the low-temperature side, we 
estimate the mass-gap ratio 1712/ mi for each N; see text for detail. In Fig. O the mass-gap ratio 
is plotted against 1/L 2 . 
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FIG. 3: The mass-gap ratio milm\ is plotted for 1/L 2 . The symbols, + and x, denote the data 
for H ^ (Fig. [2]) and H = (Fig. H|), respectively; see text for details. The least-squares fit to 
these data for the system sizes N > 7 yields m^Jmi = 1.835(4) and 1.841(4), respectively, in the 
thermodynamic limit. 




FIG. 4: Finite-size-scaling plot, (Jnn — ^viv)-^ 1 ^~ m 2/ m i> f° r ^ ne d = 3 Ising model ([2]) without 
magnetic field H = is shown. The scaling parameters are the same as those of Fig. [2j At 
Jnn — J^N = 0) we obtain the mass-gap ratio m2/mx for each system size N; see Fig. [3j 
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FIG. 5: Finite-size-scaling plot, (Jnn — J^j\f)I j ^^ u -ni2 /mi, for the d = 3 Ising model without 
magnetic field H = is shown. Here, we tentatively turned off the extended interactions Jnnn = 
Ju = 0. The critical point is set to J*nn = 0.22165455 26J], and the other scaling parameters are 
the same as those of Fig. [2J We see that the data scatter as compared to those of Fig. 01 where 
we incorporated the extended interactions. 
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